import numpy as np

a = np.random.randn(10,10)
b = np.random.randn(10)

a @ b

array([-1.06610598, -3.42056277, -4.64082216, -2.80315044,  3.59796177,
        2.13703234, -1.96828452, -0.30080257,  0.1542152 , -0.56058065])

b.T @ a.T

array([-1.06610598, -3.42056277, -4.64082216, -2.80315044,  3.59796177,
        2.13703234, -1.96828452, -0.30080257,  0.1542152 , -0.56058065])

np.eye(10)

array([[1., 0., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 1., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 1., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 1., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 1., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 1., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 1., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 1., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 1., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 0., 1.]])

np.ones((10,10))

array([[1., 1., 1., 1., 1., 1., 1., 1., 1., 1.],
       [1., 1., 1., 1., 1., 1., 1., 1., 1., 1.],
       [1., 1., 1., 1., 1., 1., 1., 1., 1., 1.],
       [1., 1., 1., 1., 1., 1., 1., 1., 1., 1.],
       [1., 1., 1., 1., 1., 1., 1., 1., 1., 1.],
       [1., 1., 1., 1., 1., 1., 1., 1., 1., 1.],
       [1., 1., 1., 1., 1., 1., 1., 1., 1., 1.],
       [1., 1., 1., 1., 1., 1., 1., 1., 1., 1.],
       [1., 1., 1., 1., 1., 1., 1., 1., 1., 1.],
       [1., 1., 1., 1., 1., 1., 1., 1., 1., 1.]])

np.zeros((10,10))

array([[0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.]])

np.linalg.norm(np.linalg.inv(a) @ a - np.eye(10))

2.570038214105557e-15

A = a

\[ Ax = b \]

x = np.linalg.lstsq(A,b)[0]

/tmp/ipykernel_2637623/3825200175.py:1: FutureWarning: `rcond` parameter will change to the default of machine precision times ``max(M, N)`` where M and N are the input matrix dimensions.
To use the future default and silence this warning we advise to pass `rcond=None`, to keep using the old, explicitly pass `rcond=-1`.
  x = np.linalg.lstsq(A,b)[0]

A @ x - b

array([-2.22044605e-15,  3.35842465e-15,  5.55111512e-15, -5.38458167e-15,
        1.77635684e-15,  2.99760217e-15,  2.66453526e-15, -2.38697950e-15,
        6.21724894e-15,  4.66293670e-15])

\[ A = U \Sigma V^T \]

\(A\) is \(m\times n\), \(U\) is \(m\times m\) \(\Sigma\) \(m \times n\) \(V\) is \(n \times n\)

\(U\) and \(V\) are orthonormal matrices, which means:

\(U U^T = I\), \(V V^T = I\)

\(\Sigma\) is 0 except for the diagonals which are the singular values of \(A\)

A = np.random.randn(100,10)

svd_A = np.linalg.svd(A)

svd_A[1]

array([13.1081352 , 11.81959892, 11.15343028, 10.83354228, 10.23738251,
        9.27330655,  9.07488313,  8.12095293,  7.844576  ,  6.77092389])

U = svd_A[0]
V = svd_A[2].T

np.linalg.norm(U @ U.T -  np.eye(100))

4.986008699316736e-15

np.linalg.norm(V @ V.T -  np.eye(10))

2.4656943504538092e-15

\[ A^+ = (A^T A)^{-1} A^T \]

\[A = U \Sigma V^T\]

\[ A^+ = (V \Sigma^T U^T U \Sigma V^T)^{-1} V \Sigma^T U^T \]

\[ A^+ = (V \Sigma^T\Sigma V^T)^{-1} V \Sigma^T U^T \]

\[ A^+ = V (\Sigma^T\Sigma)^{-1} V^T V \Sigma^T U^T \]

\[ A^+ = V (\Sigma^T\Sigma)^{-1} \Sigma^T U^T \]

\[ A^+ = V\Sigma^+ U^T \]

\[ x = V\Sigma^+ U^T b \]

\[ \Sigma^+ U^T b = V^T x \]

\[ (\Sigma^T\Sigma)^{-1}\Sigma^T U^T b = V^T x \]

\(y= V^Tx\)

\[ (\Sigma^T\Sigma)^{-1}\Sigma^T U^T b = y \]

x = V @ np.linalg.solve(S.T @ S, U.T @ b)